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ABSTRACT 

The thermal state of the intergalactic medium (IGM) provides an indirect probe 
of both the H i and He ii reionisation epochs. Current constraints on the IGM temper- 
ature from the Lya forest are restricted to the redshift range 2 < z < 4.5, limiting the 
ability to probe the thermal memory of Hi reionisation toward higher redshift. In this 
work, we present the first direct measurement of the IGM temperature around a z = 6 
quasar by analysing the Doppler widths of Lya absorption lines in the proximity zone 
of SDSS J0818+1722. We use a high resolution {R = 40 000) Keck/HIRES spectrum 
in combination with detailed numerical modelling to obtain the temperature at mean 
density, Tq = 23 600±igSj[; K (± 

9300 68 (95) per cent confidence assuming a prior 

probability 13 500 K < Tq < 38 500 K following H i and He ii reionisation. This enables 
us to place an upper limit on the redshift of Hi reionisation, zh, within 33 comoving 
Mpc of SDSS J0818+1722. If the quasar reionises the Hen in its vicinity, then in the 
limit of instantaneous reionisation we infer zjj < 9.0 (11.0) at 68 (95) per cent con- 
fidence assuming photoheating is the dominant heat source and that H i reionisation 
is driven by ionising sources with soft spectra, typical of population II stars. If the 
Hi and Hen in the IGM around SDSS J0818-I-1722 are instead reionised simultane- 
ously by a population of massive metal-free stars, characterised by very hard ionising 
spectra, we obtain a tighter upper limit of zh < 8.4 (9.4). Initiating reionisation at 
higher redshifts produces temperatures which are too low with respect to our con- 
straint unless the H i ionising sources or the quasar itself have spectra significantly 
harder than typically assumed. 
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1 INTRODUCTION 

The discovery of 2 ~ 6 quasars within the last decade has led 
to several important advances in our understanding of the 
high redshift Universe. One area where the impact of this 
work has been especially significant is the study of the high 
redshift intergalactic medium (IGM) with quasar absorption 
lines. The Lya forest in particular provides a valuable probe 
of the ionisation state of the IGM (Fan et al. 2002, 2006; 
Songaila 2004; Becker et al. 2007). However, the increasing 
opacity of the IGM to Lya photons, culminating in the ap- 
pearance of the Gunn & Peterson (1965) trough at z ~ 6 
(Becker et al. 2001; White et al. 2003), ultimately limits the 
utility of the Lya forest at the highest observable redshifts. 
This limitation has led many authors to analyse the small 
regions which exhibit transmission through the Lya forest 
even at z > 6. These highly ionised proximity zones lie be- 



tween the quasar Lya emission line and the red-most edge 
of the Gunn-Peterson trough, and are due to the enhanced 
ionisation of hydrogen close to the quasar. 

Previous analyses have focused on using these regions 
to examine the ionisation state of the IGM with low to mod- 
erate resolution spectra [R ~ 3000 — 6000). These studies 
have typically modelled the extent and/or shape of the ob- 
served transmission to obtain constraints on the IGM ion- 
isation state at z ~ 6 {e.g. Fan et al. 2006; Mesinger & 
Haiman 2007; Alvarez & Abel 2007; Bolton & Haehnelt 
2007a; Wyithe et al. 2008; MaseUi et al. 2009). However, the 
proximity zones can also be used to probe the thermal state 
of the IGM at high redshift. High resolution {R ~ 40 000) 
spectra resolve the thermal broadening kernel, enabling the 
widths of Lya absorption lines, and hence the temperature 
of the gas in the proximity zones, to be directly measured 
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(Becker et al. 2005). The line widths will be sensitive to any 
photo-heating caused by the quasar itself (Miralda-Escude 
& Rees 1994; Bolton & Haehnelt 2007b). Furthermore, since 
the cooling timescale in the low density IGM is long, infor- 
mation on the rcionisation history prior to any quasar activ- 
ity will also be encoded in these absorption lines (Haehnelt & 
Steinmetz 1998; Theuns et al. 2002a; Hui & Haiman 2003). 

Earlier measurements of the IGM thermal state using 
the Lya forest at 2 < z < 4.5 have indeed yielded valu- 
able insights into the epoch of He ii rcionisation at z ~ 3 
(Schaye et al. 2000; Ricotti et al. 2000; McDonald et al. 
2001; Zaldarriaga 2002; Lidz et al. 2009). However, given 
enough time following rcionisation, the IGM temperature 
will eventually reach an asymptotic value which depends 
only ou the spectral shape of the ionising background. In- 
formation on the earlier IGM thermal evolution, and hence 
the Hi rcionisation history at high redshift, will then be ef- 
fectively erased (Theuns ct al. 2002a; Hui & Haiman 2003). 
Consequently, it is desirable to study the thermal state of the 
IGM as close as possible to the Hi rcionisation epoch, where 
the temperature retains a more recent memory of H i photo- 
heating (Cen et al. 2009; Furlanetto & Oh 2009). The hne 
widths in proximity zones at z ~ 6 thus contain valuable 
information on the thermal history of the IGM at z>6 as 
well as the impact of any photo-heating by the quasar itself 
on its environment (Bolton & Haehnelt 2007b; Lidz et al. 
2007). 

In this paper we perform an analysis of Lya absorp- 
tion line widths in the proximity zone of the z = 6 
quasar SDSS J0818-I-1722. We compare high resolution 
Keck/HIRES data with detailed synthetic spectra con- 
structed using hydrodynamical simulations and line-of-sight 
radiative transfer to obtain the first direct constraint on the 
IGM temperature around a quasar at 2; ~ 6. As such, this 
work represents a first step toward developing a procedure 
which may be applied to larger data sets at high redshift 
in the future. We begin in section 2 by introducing the ob- 
servational data and numerical simulations wc use for our 
analysis. We discuss our methodology in section 3, and ex- 
amine the systematic uncertainties which may impact on 
our results in section 4. We present our temperature mea- 
surement in section 5. In section 6 we consider the implica- 
tions of our results for the IGM reionisation history around 
SDSS J0818+1722 before finally concluding in section 7. All 
distances are expressed in comoving units unless otherwise 
stated. 



2 DATA AND NUMERICAL MODELLING 
2.1 Observational data 

Our temperature measurements are based on a 
high-resolution spectrum of the 2 — 6.00 quasar 
SDSS J0818-hl722 (Fan ct al. 2006). Keck HIRES 
data were taken in February 2006 using the upgraded 
detector. We employed a 0.86"slit, which produces a 
resolution oi R = 40 000 (FWHM = 6.7 km s"^). The 
total integration time was 7.5 hours. The raw data were 
processed using a custom set of IDL routines that includes 
optimal sky subtraction (Kelson 2003). In order to achieve 
the highest possible signal-to-noise (S/N) ratio, the final 



one-dimensional spectrum was optimally extracted from 
all exposures simultaneously after relative flux calibrations 
were applied to remove the blaze functions from individual 
orders. This allowed us to efficiently reject cosmic rays and 
other bad pixels while preserving as many counts from the 
object as possible. The final S/N near the Lya emission 
line was ~ 15 per 2.1 km s~^ binned pixel. 

Continuum fitting in the proximity zone region of a 
z ^ 6 quasar is challenging due to the strong absorp- 
tion present in the blue side of the Lya emission line. We 
first divided the combined spectrum by a power-law with 
Fi, oc i/~ ■ , normalised near 1280 (1 + z) A. The value of 
the power-law slope only weakly affects the overall contin- 
uum level in the proximity zone. Next wc fit a slowly varying 
spline to the emission line. The blue side of the line profile 
was made to roughly mirror the red side, which could be 
drawn though the unabsorbed continuum. This produced a 
generally Gaussian shape, which wc joined smoothly onto 
the underlying power-law. The continuum was set ~ 10 
per cent above the tops of the transmission peaks within 
~ 2000 km s~^ of the quasar redshift, after which the peats 
were allowed to fall further below the continuum. Although 
the continuum level is obscured by the high levels of absorp- 
tion, we estimate that our fit is within 20 per cent of the 
correct value over the wavelength range of interest. 

2.2 Hydrodynamical simulations 

The synthetic proximity zone spectra used in this study are 
constructed using line-of-sight density, velocity and temper- 
ature fields drawn from cosmological hydrodynamical simu- 
lations. The simulations were performed using a customised 
version of the parallel Tree-SPH code GADGET-3, which is 
an updated version of the publicly available code GADGET- 
2 (Springel 2005). We use twelve different hydrodynamical 
simulations in this work, summarised in Table 1. Nine simu- 
lations (A-H and W) were performed in 10h~^ Mpc periodic 
boxes containing 2 x 512^ gas and dark matter particles. The 
simulations are specifically designed to resolve the Lya for- 
est at high redshift (Bolton & Becker 2009). However, they 
employ a relatively small box size to achieve the required 
mciss resolution. In order to assess the effect of box size 
on our results, we performed two further simulations with 
identical mass resolution but different box sizes of Wh~^ 
Mpc (LI) and 40/i~^ Mpc (L3), respectively. Finally, another 
10ft~^ Mpc simulation, L2, was performed with 2 x 256^ gas 
and dark matter particles to provide an additional test of 
convergence with mass resolution. 

The simulations were all started at 2; = 99, with 
initial conditions generated using the transfer function of 
Eisonstcin & Hu (1999). The cosmological parameters are 
fim = 0.26, = 0.74, fibft^ = 0.023, h = 0.72, as = 0.80, 
Us = 0.96, consistent with recent studies of the cosmic mi- 
crowave background (Komatsu ct al. 2009; Roichardt et al. 
2009). The IGM is assumed to be of primordial composition 
with a helium fraction by mass of Y = 0.24 (Olive & Skill- 
man 2004). The gravitational softening length was set to 
1/30*'' of the mean linear interparticle spacing and - with 
the exception of model W - star formation was included 
using a simplified prescription which converts all gas par- 
ticles with overdensity A = p/(p) > 10^ and temperature 
T < 10^ K into collisionless stars. Model W instead uses the 
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Table 1. Hydrodynamical/raxliative transfer simulations used in this work. Prom left to right, the columns list the simulation identifier, 
the box size, the total particle number, the gas particle mass, the scaling factors for the UVB photo-heating rates (see main text for 
details) , the median volume weighted gas temperature at mean density. To , and the power law slope of the temperature-density relation, 
T = TqAT"!, both prior to and after Hen photo-heating by the quasar. The slope is given in terms of 7 — 1 = 2(logTo — logT-o.s), 
where T_o.5 is the median volume weighted gas temperature at log A = —0.5. The temperatures are quoted to three significant figures 
only. Note a power-law temperature-density relation following Hell photo- heating by the quasar is an approximation only (see Figure 1). 



Model L[/i-iMpc] Particles Mgas [/i^^M©] C i To,i [K] 7; - 1 To [K] 7-I 



A 


10 


2 X 512^ 


9.2 X Iff* 


0.30 


0.0 


4 300 


0.39 


13 500 


0.12 


B 


10 


2 X 5123 


9.2 X 10* 


0.80 


0.0 


8 500 


0.40 


17500 


0.18 


C 


10 


2 X 512^ 


9.2 X 10" 


1.45 


0.0 


13 000 


0.41 


21500 


0.22 


D 


10 


2 X 512^ 


9.2 X Iff* 


2.20 


0.0 


17300 


0.41 


25 700 


0.25 


E 


10 


2 X 5123 


9.2 X 10* 


3.10 


0.0 


21900 


0.41 


29 900 


0.27 


P 


10 


2 X 5123 


9.2 X 10" 


4.20 


0.0 


26800 


0.41 


34 400 


0.28 


G 


10 


2 > 


< 5123 


9.2 X Iff* 


5.30 


0.0 


31100 


0.40 


38 500 


0.29 


H 


10 


2 > 


< 5123 


9.2 X Iff* 


1.45 


-1.0 


12 800 


0.06 


21400 


0.04 


LI 


10 


2 > 


< 1283 


5.9 X 10^ 


1.45 


0.0 


13 600 


0.40 


22 300 


0.22 


L2 


10 


2 > 


< 2563 


7.3 X IQS 


1.45 


0.0 


13 200 


0.40 


22 000 


0.21 


L3 


40 


2 > 


< 5123 


5.9 X 10« 


1.45 


0.0 


13 400 


0.41 


20 500 


0.23 


W 


10 


2 X 5123 


9.2 X 10* 


1.45 


0.0 


13 000 


0.41 


21700 


0.22 



multi-phase star formation and feedback model of Springel 
& Hernquist (2003). We use this to explore the impact of 
strong galactic winds on our results. 

The gas in the simulations is assumed to be optically 
thin and in ionisation equilibrium with a spatially uniform 
ultraviolet background (UVB). The UVB corresponds to the 
galaxies and quasars emission model of Haardt & Madau 
(2001). Hydrogen is reioniscd at z = 9 and gas with A < 10 
subsequently follows a tight power-law temperature-density 
relation, T = TqA''"^, where To is the temperature of the 
IGM at mean density (Hui & Gnedin 1997; Valageas et al. 
2002). In order to explore a variety of thermal histories, we 
rescale the Haardt & Madau (2001) photo-heating rates by 
different constants in models A-H. In each simulation we as- 
sume ei = CA^e"'''"'\ where are the Haardt & Madau 
(2001) photo-heating rates for species i = [H i, He i. He 11] 
and ^ are constants listed in Table 1. The LI, L2, L3 and 
W simulations use same thermal history as model C. A den- 
sity dependent roscaling of the photo-heating rates has been 
applied to model H, giving 7 ~ 1 while maintaining a similar 
To to model C. All the other hydrodynamical simulations in 
this study have 7 ~ 1.4. Simulation outputs were obtained 
from each model at z = 6.01. The values of To and 7 in the 
hydrodynamical simulations are listed in columns 7 and 8 of 
Table 1. 



2.3 Line-of-sight radiative transfer models 

In order to correctly model photo-ionisation and heating^ 
by the quasar, we have implemented a line-of-sight radiative 



^ The ionisation and thermal state of the IGM is fully recom- 
puted in post-processing by our radiative transfer algorithm; 
we only use the hydrodynamical simulations for the initial con- 
ditions (gas density, peculiar velocity and temperature) of the 
gas. The gas overdensity and peculiar velocity field thus re- 
main static in the RT calculation, and the hydrodynamical re- 
sponse of the gas to the photo-heating induced by the quasar is 
not modelled. This effect should be small over the quasar life- 
time we consider (tq ~ 10^ yr, cf. the sound crossing timescale 



transfer (RT) scheme in the hydrodynamical simulations to 
account for the propagation of ionising radiation into the 
IGM. For this we use the multi-frequency line-of-sight rar 
diative transfer algorithm described by Bolton & Haehnelt 
(2007b) and updated in Bohon et al. (2009). 

We follow the procedure outlined in Bolton & Haehnelt 
(2007b) for constructing the initial conditions for our sim- 
ulations including radiative transfer. We use a friends-of- 
friends halo finding algorithm with a linking length of 0.2 to 
identify haloes in the hydrodynamical simulations. We then 
select the ten most massive haloes in each simulation and 
extract lines-of-sight in diflferent orientations around them. 
The exact halo mass has little impact; wo find it is more 
important to resolve the IGM along the line-of-sight rather 
than the host halo itself. We discuss this issue further in sec- 
tion 4.1. The halo lines-of-sight are then spliced with lines- 
of-sight drawn at random from the rest of the simulation, 
resulting in 30 different lines-of-sight 55h~^ Mpc in length, 
all of which start at the location of a halo. 

For each line-of-sight we compute the transfer of ionis- 
ing radiation through the IGM from a quasax which emits 
N photons per second above the Hi ionisation threshold, 



N 



dv. 



(1) 



L,j oc 



-0.5 



(2) 



with a generic broken power law spectrum given by 

(1050 < A < 1450 A), 
(A < 1050 A). 

This spectrum has an extreme UV (EUV) spectral index 
as = 1-5, consistent with radio quiet quasars at lower red- 
shift (Tolfor et al. 2002). Note, however, the exact spectrum 
and hence ionising luminosity of the quasar is rather un- 
certain at 2 = 6. Adopting a harder (softer) EUV spectral 
index will increase (decrease) the ionising luminosity and 
the amount of photo-heating around the quasar. We assume 



tsc ^ 6.7 X 10*(L/10properkpc)(T/10'* K)-V2yr and the dy- 
namical timescale tdyn ^ 4.2 x IQ^ + z)/!]-^/"^ yv). 
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A4"i45o = —27.4 for the absolute magnitude of the quasar, 
corresponding to iV = 2.74 x lO^'^s"^ for our chosen spec- 
trum. This magnitude corresponds to the M1450 obtained 
for SDSS J0818+1722 (Fan et al. 2006). Lastly, we use a 
duration of Iq — 10^ yr for the phase when the quasar is 
optically bright (Haehnelt et al. 1998; Martini 2004; Cro- 
ton 2009). This parameter is also somewhat uncertain; we 
discuss the implications of this in more detail shortly. 

In all our RT simulations, we initialise the ionisation 
state of the IGM by assuming photo-ionisation equilibrium 
with a UVB given by 



Jl/ = J-21 (—^ X 



(fHI <l^< J^Hell) 



(3) 



where J-21 is the amplitude of the UVB at the H i ionisation 
threshold in proper units of lO"^'^ ergs~^ cm~^ Hz~^ sr~^. 
We adopt a soft, stellar like spectrum and assume Hi and 
He I is already highly ionised by z < 6 {e.g. Choudhury 
& Ferrara 2006; Becker et al. 2007; Pritchard et al. 2009). 
There is no evidence to suggest the hydrogen in the IGM sur- 
rounding SDSS J0818+1722 was significantly neutral prior 
to the quasar turning on (see Bolton & Haehnelt 2007b 
for a detailed discussion of this point). The exact choice of 
UVB spectral index at frequencies below the He 11 ionisation 
threshold is unimportant, since the IGM temperature is ini- 
tialised using the hydrodynamical simulations. 

However, Hen reionisation, which must be driven by 
sources with hard ionising spectra, is thought to be delayed 
until around 2; ~ 3 when the number density of quasars 
begins to peak (Madau & Meiksin 1994; Furlanetto & Oh 
2008). Wc model this delay by truncating the UVB above the 
He II ionisation threshold. When the quasar turns on in our 
RT models, its hard ionising spectrum will therefore reionisc 
and photo-heat the He 11 in the surrounding IGM, resulting 
in the so-called "thermal proximity effect" (Miralda-Escude 
& Rees 1994; Thcuns ct al. 2002b; Meiksin et al. 2010). 
The main effect of He 11 reionisation by the quasar is to flat- 
ten the temperature density relation and boost the tem- 
perature at mean density by ~ 7000 — 9000 K within the 
proximity zone. The volume weighted temperature-density 
plane from one of our RT simulation sets (model C) is dis- 
played in Figure 1. The final (approximate) temperature- 
density relation in the RT simulations after He 11 photo- 
heating by the quasar are given in the last two columns 
of Table 1. As already discussed, the exact amount of heat- 
ing depends on the assumed EUV spectral index. Note also 
that if the quasar lifetime is significantly shorter (longer) 
than the generic value of Iq = 10^ yr we have assumed, the 
extent of the region where He 11 is photo- ionised and heated 
by the quasar will be smaller (larger). In this work we con- 
sider the IGM within 33 Mpc of SDSS J0818-hl722, a scale 
over which He 11 is fully ionised for Iq = 10^ yr for our choice 
of quasar spectrum. Longer quasar lifetimes will make little 
difference to our results. However, if Iq < 10^ yr our mod- 
els will overestimate the size of the He iii region (and hence 
Hen photo-heating) around the quasar. 

Lastly, the amplitude of the UVB at the H i ionisation 
edge in Eq. (3), J-21, is set using the density dependent 
model presented in Wyithc ct al. (2008). Briefly, this model 
accounts for the bias in the overdense region around the 
quasar host halo and is computed as a function of the proper 
time along the trajectory of a photon emitted by the quasar. 




Figure 1. Scatter plot of the volume weighted temperature- 
density plane for the 30 lines-of-sight used in model C. A random 
subset of pixels 250 km s^-' < vh < 3250 km s^-' from the quasar 
redshift are displayed for clarity of presentation. The dashed line 
shows the final (approximate) power-law temperature-density re- 
lation after Hen photo- heating by the quasar, while the dotted 
line corresponds to the initial IGM temperature-density relation 
(taken from the hydrodynamical simulation) before the quasar 
turned on. There is significant scatter in the relationship between 
temperature and density following quasar activity. The shaded 
region corresponds to 95 per cent of all pixels with normalised 
flux 0.05 < F < 0.95 in the corresponding proximity zone spec- 
tra. This gives an indication of the densities which transmission 
in the proximity zone is sensitive to (see the discussion in sec- 
tion 4.2). 



This results in an enhanced contribution to the UVB near 
the quasar host halo relative to the mean IGM. The lat- 
ter is calibrated to match measurements of the Hi photo- 
ionisation rate derived from the observed Lya forest opacity 
at 4 < 2: < 6 (Bolton & Haehnelt 2007c). A detailed de- 
scription of the model and its application to high redshift 
Lya forest spectra may be found in Wyithe et al. (2008). 
Note that the UVB amplitude is important for the observed 
size of a proximity zone (Wyithe et al. 2008), but the pre- 
cise value of this parameter is not important for the gas 
temperatures around the quasar we explore in this work. 



2.4 Spectra construction 

We follow standard procedure to construct synthetic prox- 
imity zone spectra from the output of each of our RT mod- 
els {e.g Thcuns ct al. 1998). The spectra are then con- 
volved with a Gaussian instrument profile with FWHM = 
6.7kms~^, resampled onto pixels of width 2.1kms~^ and 
Gaussian distributed noise with S/N = 12 is added. 
These values are chosen to match the observed spectrum of 
SDSS J0818+1722. To account for uncertainties in the con- 
tinuum placement on the observed data, we renormalise the 
synthetic data by the highest flux in 1000 km s~^ segments 
along each line-of-sight. Although this will not account for 
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Figure 2. Upper four panels: The topmost panel shows an example synthetic spectrum at z = 6 drawn from model C. The quasar is 
situated at t'u = km s^^ on the right hand side of each panel. The vertical dashed lines indicate the velocity range over which Voigt profiles 
arc fitted to the data, 250 kms^^ < vjj < 3250 kms^-"^ (2.6Mpc < R < 33.3Mpc). The vertical arrows mark the redshift of the Voigt 
profile fits made with VPFIT, while the fits themselves are shown by the dotted curves. The corresponding gas temperature, overdensity 
and fractional abundances of Hi and Hen along the line-of-sight are shown in the subsequent three panels. Lower panel: Keck/HIRES 
spectrum of SDSS J0818+1722 at z = 6.00. 



© 0000 RAS, MNRAS 000, 000-000 



6 J.S. Bolton et al. 




20 40 60 80 20 40 60 80 100 

b [km s"'] b [km s"'] 



Figure 3. Left: The Dopplcr parameter CPDF for four of our RT models. The simulated CPDFs correspond to the data from all 
30 lines-of-sight in each model. The CPDF is clearly sensitive to the gas temperature in the quasar proximity zone, with the CPDF 
shifting to higher b-valucs for increasing temperature. The median Doppler parameters for the simulated models displayed are femed = 
(24.5, 30.0, 33.2, 38.4)kms"l. Right: The observed distribution for SDSS J0818+1722 (thick solid curve), which has ftmcd = 31.0kms-\ 
in comparison to the CPDFs for the individual lines-of-sight in model C. This gives an indication of the scatter in the CPDF from one 
line-of-sight to the next. 



continuum uncertainties perfectly, it will approximate the 
way the continuum is fitted to the real data and should 
reduce a potential systematic bias. We consider this issue 
further in section 4.2. Lastly, we rescale the optical depths 
in each pixel of the synthetic spectra by a constant to match 
the moan flux of the observed data, {F)obs = .^obs/.fcont = 
0.486, within the range 250 km s"^ < vh < 3250 km s"^ 
(2.6 Mpc < R < 33.3 Mpc). An example spectrum drawn 
from model C is displayed in Figure 2 along with the ob- 
served Une-of-sight we use in this work, SDSS J0818-I-1722. 



3 ANALYSIS PROCEDURE 
3.1 Voigt profile fitting 

In this work wc use the cunmlative probability distribution 
function (CPDF) of the Doppler widths of Lya absorption 
lines, 6, in the proximity zone as our probe of the IGM tem- 
perature. The CPDF has the advantage of fully using the 
limited data available (our analysis of SDSS J0818-I-1722 
yields only 25 Dopplcr widths) and it avoids binning and the 
associated loss of information. Although the absorption lines 
will not all be purely thermally broadened (Theuns et al. 
2000), all absorption lines will nevertheless be smoothed 
on a scale associated with the thermal broadening kernel. 
The entire CPDF is therefore sensitive to the IGM temper- 
ature. However, the lack of a direct relationship between the 
Doppler parameters and temperature means detailed syn- 
thetic spectra are crucial for calibrating this statistic. 

We obtain the CPDF for both the observed and syn- 
thetic data via Voigt profile fitting. We perform the anal- 
ysis using an automated version of the Voigt profile fitting 



package VPFIT.^ We choose to include lines in the range 
250 km s"^ < VH < 3250 km s"^ (2.6 Mpc < R < 33.3 Mpc) 
from the quasar redshift only. The lower limit is chosen 
to avoid edge effects in the synthetic data, while beyond 
3250 km s"^ it becomes impossible to fit Voigt profiles reli- 
ably due to the increasing opacity of the IGM. All lines with 
relative errors in excess of 50 per cent are rejected from our 
final analysis (Schaye et al. 1999). This removes most of the 
very narrow lines which tend to be in blends. Wc also dis- 
card all lines with Dopplcr parameters b > 100 km s^^ and 
column densities log A^hi > 17. Broad shallow lines with 
b > 100 km s"'^ are almost always added near the contin- 
uum level to improve the overall fit, while strong lines with 
damping wings tend to be added near the edge of the prox- 
imity zone where there is very little transmission. These cuts 
remove 32 and 25 per cent of identified Voigt profiles in the 
observed and simulated data sets, respectively. Lastly, we 
note that some ambiguous fits will still remain in our analy- 
sis. Ifowever, since we treat the observed and synthetic spec- 
tra in exactly the same way we expect any potential bias to 
be small. 

The sensitivity of the resulting Doppler width CPDF 
to the gas temperature is demonstrated in the left panel of 
Figure 3, where the CPDFs for four of our RT simulation 
sets are displayed. The CPDFs correspond to the data from 
all 30 lines-of-sight in each RT simulation set. It is clear the 
entire CPDF is sensitive to the gas temperature in the prox- 
imity zone. Increasing the temperature of the IGM shifts 
the CPDF towards higher b-parameter values, increasing the 
median Doppler parameter. The median Doppler parameters 
for the models displayed, which spam a range of 13 500 K < 
To < 38500K, are b^ed = (24.5, 30.0, 33.2, 38.4) kms'^ 

^ http: / /www. est . cam. ac. uk /~rfc/ vpfit .html 
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The observed distribution, which has 6med ~ 31.0 km s^^, 
is displayed as the thick solid curve in the right panel of 
Figure 3 along with the CPDFs for the individual lines-of- 
sight in model C. This gives an indication of the scatter in 
the CPDF from one linc-of-sight to the next for the same 
underlying model. A two distribution Kolmogorov-Smirnov 
(KS) test (Press et al. 1992) using the average (all 30 lines- 
of-sight) and observed CPDFs yields Pks = 0.06, 0.94, 0.21 
and 4x 10"'' for models A, C, E and G respectively. This does 
not rule out the null hypothesis that the observed CPDF is 
drawn from the same distribution as models A, C and E. 

3.2 Monte-Carlo analysis 

We now turn to describing the methodology which lies at 

the heart of our analysis procedure. We use our synthetic 
spectra to construct Monte Carlo realisations of the Dopplcr 
parameter CPDF for a range of models with different IGM 
temperatures. To quantify the amount of scatter in the sim- 
ulated CPDFs, we use a "D-statistic" that is very similar 
to the parameter used in a KS test. For each linc-of-sight 
in a given simulation set, the D-statistic is the maximum 
difference between the Doppler parameter CPDF for that 
line-of-sight and the CPDF for all 30 lines-of-sight, such that 

A = max|P(< 6)i - P(< 6)aii|, i = 1...30, (4) 

where we preserve the sign of the difference. The D-statistic 
CPDF for a model with known temperature To, P{< D\To), 
can then be constructed. If the D-statistic for an observed 
linc-of-sight. Dobs, is known, the D-statistic CPDF may then 
be used in combination with Bayes theorem to infer a confi- 
dence interval for the observed Tq. The cumulative probabil- 
ity for the observed temperature at mean density To given 
Dobs is: 

P« TolAbs) = £° dP{<U\Dob.) ^^, 

An application of Bayes theorem, p(To|-Dobs) oc 
p{Dohs\To)p{To), then leads to 

where K is & constant which normalises the total probability 
to unity and p{To) is the prior on To. In order to obtain our 
constraint on To, we evaluate the derivative of the D-statistic 
CPDF derived from our simulations at Dobs in Eq. (6). We 
adopt a flat prior, p(To), over the range of temperatures ex- 
plored in the simulations. These are set by the combined 
effect of the initial temperature of the hydrodynamical sim- 
ulations and the subsequent heat input from Hen photo- 
heating by the quasar, such that 13 500 K < To < 38 500 K. 
This is intended to represent a reasonable range for the 
IGM temperature following Hi and Hcii rcionisation, and 
is consistent with temperatures predicted in simulations of 
Hen rcionisation (McQuinn et al. 2009). 

The main drawback of our method is the computational 
expense of constructing and analysing synthetic spectra. We 
have seven hydrodynamical simulations with different values 
of To, so we may only evaluate ''"'''''"^o'"'' ^""' at seven discrete 
points. We use a cubic spline to interpolate between these 
to obtain a continuous distribution for integration. With the 



current method, finer sampling of the distribution would re- 
quire running and analysing additional simulations. A less 
expensive approach would be to use a single hydrodynamical 
simulation and simply impose a range of initial temperature- 
density relations on the IGM. However, this would decouple 
the initial gas temperature from the gas hydrodynamics and 
incorrectly model the effect of pressure (Jeans) smoothing 
on the gas distribution {e.g. Pawlik et al. 2009; Peeples et al. 
2009). Furthermore, there are only 30 line-of-sight in each 
simulation set, which makes the derivative of the D-statistic 
CPDF somewhat noisy. This could be remedied by analysing 
more lines-of-sight, although the Voigt profile analysis pro- 
cedure is time consuming. For this first attempt we have 
instead elected to perform a easily manageable number of 
simulations. We therefore smooth the D-statistic CPDF with 
a Gaussian filter of width a = 0.025 before computing its 
derivative. We find this to be the optimum smoothing width, 
and we have verified that alternative choices of a = 0.01 and 
0.05 do not significantly alter our final results. 

On the other hand, the major advantage of our ap- 
proach is that it avoids binning the data and is non- 
parametric. However, our method implicitly assumes that 
the synthetic spectra accurately represent the observed data. 
If our synthetic models are significantly in error, any tem- 
perature constraints will be unreliable. For this reason, we 
now turn to consider numerical convergence and some of the 
possible systematic errors which may affect our analysis. 



4 CONVERGENCE AND SYSTEMATICS 
4.1 Numerical convergence 

We first consider the effect of mass resolution and box size 
on the simulated Doppler width CPDF in Figure 4. The 
left panel displays the impact of mass resolution on the 
CPDF, while the right panel indicates the effect of box size. 
The CPDF is well converged for our fiducial mass resolu- 
tion (9.2 X 10* Mq), although note our lowest resolution 
model (LI, 5.9 x W^h^^ M©) does not correctly resolve the 
line widths. It is also encouraging that the box size, and 
hence the amount of large scale power in our simulations, 
has little effect on the CPDF. The maximum halo mass in 
the 40h~^ Mpc box is 4.33x 10"/i"^ M©, which is an order of 
magnitude greater than the largest halo in our fiducial simu- 
lation volume of lOft"^ Mpc (6.13 x 10^°/i"^ Mq). However, 
this has little effect on the CPDF because we arc measur- 
ing Doppler parameters in the ambient IGM rather than in 
the host environment of the quasar. Note, however, that the 
maximum halo meiss we consider is still smaller than the 
expected host masses of ^ 10^^ — 10"'^'' Mq for bright z ~ 6 
quasars (Walter et al. 2004; Volonteri & Recs 2006; Li et al. 
2007). We nevertheless find it is much more important to 
resolve absorption line widths correctly in our models using 
a high resolution simulation rather than using a larger vol- 
ume. We conclude that our simulations arc well converged 
with box size and mass resolution for the purpose of this 
study; any differences in the CPDF are small compared to 
the contrast between models with different temperatures. 
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Figure 4. Numerical convergence tests for tlie simulated Doppler width CPDF. The CPDFs correspond to the data from all 30 lines-of-sight 
in each RT simulation set. Left: Test of the hydrodynamical simulation mass resolution. Although the CPDF from the lowest resolution 
LI model has not converged, the distribution is well captured for our fiducial mass resolution of 9.2 X IO-^/i-^Mq (model C) Right: Test 
of the simulation box size. The CPDF is well converged, indicating our choice of a 10h~^ Mpc box is sufficient. 



4.2 Systematic uncertainties 

4-2.1 The mean flux 

The effect of four different systematic uncertainties on the 
simulated Doppler parameter CPDF are displayed in Fig- 
ure 5. The first we consider is the mean transmitted flux, 
{F) = (e""^'), where n is the optical depth in pixel i. As 
discussed earlier, when constructing our synthetic spectra 
we follow a standard procedure and rescale the Hi photo- 
ionisation rate in post-processing by a fixed constant to pro- 
duce spectra which match the (F) of the observed data, 
{F) = 0.486. It is therefore important to assess how a dif- 
ferent value of (F) may influence the Doppler parameter 
CPDF. 

The effect of the mean flux in the region of the spec- 
tra where lines are fitted, 250 km s"'^ < vh < 3250 km s~^, 
is shown in the upper left panel of Figure 5. We perform 
a Voigt profile analysis on two further models which are 
identical to model C, except from being roscaled in post- 
processing to have (F) — 0.386 and (F) = 0.586. It is clear 
that the mean fiux has little effect on the Doppler parameter 
CPDF, changing the median Doppler parameter by less than 
0.5 kms~^ for the examples shown. This is expected theoret- 
ically; as pointed out by Hui & Rutlcdgc (1999), the ionising 
background does not play a role in setting the shape of an 
absorption line - and hence its Doppler parameter - axound 
the peak optical depth. Note, however, that the mean fiux 
will indeed be important for the line column densities, de- 
creasing A^Hi at fixed density when the photo-ionisation rate 
is raised. We find that uncertainties in the mean fiux should 
not be important for this work. 



4-2.2 Galactic winds 

We reach a similar conclusion regarding the impact of galac- 
tic winds on the CPDF in the upper right panel of Fig- 



ure 5. We have constructed spectra from model W. This uses 
the multi-phase star formation and winds model of Springel 
& Hernquist (2003), but in all other respects is identical 
to model C. The wind velocity is Vy, = 448kms~^, which 
is an extreme model in the sense that it assumes all en- 
ergy produced by supcrnovac is converted into kinetic en- 
ergy The differences in the CPDF are again small, with 
&med = 30.0 kms~^ and 30.3 kms~^ for the C and W mod- 
els respectively. 



4-2.3 The temperature- density relation 

We do not attempt to measure the density dependence of 
the IGM temperature in this work; this is difficult to achieve 
for a single quasar spectrum at z = 6. Instead we focus on 
constraining the temperature at mean density. To- The as- 
sumed power-law slope of the temperature-density relation, 
7, is therefore a systematic uncertainty in our analysis. Note 
again that a single power-law temperature-density relation, 
T = TqA"'^^ , docs not provide an accurate description of 
the IGM thermal state in our simulations (e.g. Figure 1, 
and see also Bolton et al. 2004; Trac et al. 2008). Never- 
theless, for the present we find it convenient to characterise 
the relationship between temperature and density as scatter 
around a power-law. 

In the lower left panel of Figure 5 we demonstrate the 
effect of 7 on the Doppler parameter CPDF. Model C and 
model H have 7 ~ 1.2 and 7 ~ 1.0 respectively, with almost 
identical temperatures at mean density. Interestingly, there 
is very little difference between the CPDFs, with femed = 
30.0 kms~^ and 29.3kms~'^. This is the result one would 
expect if the CPDF is sensitive to gas densities around the 
cosmic mean at z = 6. We may test if this is the case by 
evaluating the range of optical depth weighted overdensities. 
At- = ^ AiTi/ ^ n (e.p. Schaye et al. 1999), over the range 
250 km s"^ < Vh < 3250 km s in our simulated spectra. In 
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Figure 5. Systematic uncertainties and their eflFect on the Doppler parameter CPDF. Upper left: The effect of the mean flux measured 
over 250kms~^ < vh < 3250 kms"^. The fiducial value in this work is {F) = 0.486 (solid curve), corresponding to value measured from 
the spectrum of SDSS J0818+1722. Two further models with (F) = 0.386 (dashed curve) and (F) = 0.586 (dotted curve) are displayed for 
comparison. Upper right: The impact of feedback in the form of strong galactic winds, Vw = 448kms~^, implemented using the multi-phase 
star formation model of Springel & Hernquist (2003) . Model C (solid curve) has no winds while model W (dotted line) is identical to model 
C in all respects aside from the inclusion of galactic winds. Lower left: The effect of the temperature-density relation slope, 7. Model C 
(solid curve) has 7 ~ 1.2 while model H (dotted curve) has 7 ~ 1 following Hell photo-heating by the quasar. Lower right: The effect of 
uncertainty in the continuum placement. The dotted curve shows the CPDF obtained before the standard continuum correction is applied 
to the synthetic data (see main text for details), while the solid curve is obtained from the data after correction. Lowering the continuum 
(it is never raised in the correction) always reduces the median Doppler parameter. 



Figure 1 the shaded region corresponds to the range of 
including 95 per cent of the pixels with 0.05 < F < 0.95 in 
the proximity zone: 0.4 < < 4.4, with a median value 
of Ar = 0.95. The transmission in the proximity zone thus 
predominantly probes the IGM around mean density, and it 
is clear why different values of 7 have only a small impact on 
the CPDF. Fiiture studies at z ~ 6 could utilise higher order 
Lyman series absorption sensitive to higher gas densities to 
obtain constraints on 7 (Dijkstra et al. 2004; F^rlanetto & 
Oh 2009). For now, based on these results we shall assume 
the impact of this parameter on our analysis is small. 



4.2.4 Continuum placement 

We now consider the effect of the uncertainty in the intrin- 
sic emission of the quasar, Jcont, where Jobs ~ Jconte"^. At 
z = 2 — 3, the typical uncertainty in the continuum place- 
ment in high resolution Lya forest data is estimated to be 
1-4 per cent (Tytlcr ct al. 2004; Kirkman et al. 2005; Kim 
et al. 2007; Faucher-Gigucrc et al. 2008). However, towards 
higher redshift, the increasing opacity of the IGM means the 
transmitted flux recovers to the unabsorbed level less and 
less frequently. Identifying the continuum level accurately 
becomes very challenging, with an estimated uncertainty of 
up to 20 per cent even for the high-quality spectrum con- 
sidered here. Continuum estimation therefore presents a sig- 
nificant systematic uncertainty in our analysis. 
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Figure 6. Test of our analysis procedure using two random lincs-of-sight drawn from our simulated data set. Tiie lines-of-sigiit in the left 
and right panels are drawn from models B and F, respectively. The solid curves show the cumulative distribution for To, P(< TolDobs)! 
computed using Eq. (6). The filled circles show the "true" temperature of the models at mean density, while the shaded regions indicate 
the 68 per cent confidence interval around the recovered median. The numerical values in the top left of each panel correspond to the 
recovered median temperature and the 68 (95) per cent confidence intervals. The large uncertainties are dominated by the statistical error 
associated with measuring the temperature along a single line-of-sight. 



As it is difficult to identify where the continuum lies 
on the observed data at z ~ 6, our approach in this work 
is instead to treat the simulated spectra in the same way 
as the observational data. The observational data were nor- 
malised using a slowly varying spline fit, adjusted by hand 
to roughly follow the profile of the transmission peaks in 
the proximity zone. The continuum correction we apply to 
our data, described in section 2.4, is designed to mimic this. 
To recap, wc rcnormalisc the synthetic data by the highest 
fiux in ~ 1000 km s^^ segments along each line-of-sight. Al- 
though this means wc may not correctly identify the "true" 
continuum, in principle this procedure should at least reduce 
any bias (see Ranch et al. 1997; McDonald et al. 2000 for 
similar approaches to this problem at lower rcdshift). 

The lower right panel of Figure 5 displays the effect 
of the continuum correction on the CPDF, which always 
acts to lower the continuum on the synthetic spectra. Al- 
though this correction has little effect near the redshift of 
the quasar - the IGM here is very highly ionised and the flux 
regularly reaches the unabsorbed level in our simulations - 
it becomes larger towards the blue-most edge of the region 
we fit lines to; the Lya transmission rarely reaches F = 1 
beyond ~ 2000 kms^^. On average the maximum correction 
is around 10 per cent, but it can be up to a 25 jjcr cent for 
some lines-of-sight at the edge of the proximity zone. The 
correction shifts the CPDF to lower Doppler widths, effec- 
tively narrowing the absorption lines as the continuum is 
lowered. For reference, the median Doppler parameters are 
6med = 30.0 kms~^ (with correction) and 31.6 kms^^ (no 
correction). Consequently, if the continuum on the observed 
data has been placed high (the correction always lowers the 
continuum) relative to the corrected synthetic data, wc will 
infer gas temperatures which are ~ 2000 K too high. This is 
clearly an important systematic to consider when the statis- 
tical error bars approach the level of a few thousand degrees. 



4.2.5 Metals 



Lastly, wc consider the impact of mis-idcntificd metal lines 
on the Doppler parameter CPDF. A full treatment would 
require modelling the distribution and abundance of metals 
in the IGM at a = 6, which is beyond the scope of this work. 
However, assuming metal lines are typically narrower than 
the Lya lines, wc may estimate the effect of any metal con- 
tamination as follows. In our spectrum of SDSS J0818-I-1722, 
there are roughly ~ 2 metal absorption lines per 3000 km s~^ 
interval redward of the Lya emission line. Over the same in- 
terval in the proximity zone, by comparison, we include 25 
lines in our Doppler parameter CPDF. We may therefore 
estimate the effect of possible metal contamination by ex- 
cluding the two narrowest lines identified in the proximity 
zone of SDSS J0818-I-1722 from our analysis, 6 = 8.6 kms"^ 
and h = 10.0 kms^^ (c/. the narrowest line in model C, 
b = 7.9 km s^^). This changes the median Doppler param- 
eter for SDSS J0818-I-1722 from 6^ed = 31.0 km s"^ to 
&med = 31.3 km s~^. We estimate that metal contamination 
may then at most systematically decrease any constraint on 
the IGM temperature in the proximity zone by ~ 2000 K. As 
we shall see, the systematic uncertainties discussed in this 
section are small compared to the statistical uncertainty for 
an individual line-of-sight. Consequently, we do not consider 
these in our final analysis. However, we note that future 
analyses with smaller statistical errors will require a more 
detailed treatment of continuum placement and metal con- 
tamination. 
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5 RESULTS 

5.1 Test of the analysis procedure 

We arc now ready to progress to the main result of this 
paper. However, before proceeding further we briefly test our 
methodology. We apply the analysis procedure described in 
section 3 to two synthetic lines-of-sight for which we already 
know the "true'" temperature. In this way we may check if 
our procedure correctly recovers the temperature of the IGM 
around a quasar. 

The test of our methodology is displayed in Figure 6. 
Each panel shows the constraint on To for a synthetic line- 
of-sight in the form of a cumulative distribution for To, P{< 
To j -Dobs), computed using Eq. (6). The lines-of-sight are 
drawn from the B and F models which have Tq''"'' = 17 500 K 
and Tq^"" = 34 400 K, respectively. These temperatures are 
represented by the filled circles in Figure 6. The shaded 
range in both panels display the 68 per cent confidence in- 
tervals around the inferred median value of To. Due to the 
limited statistical power of a single line-of-sight the uncer- 
tainties on the constraints are large, Tq = 18 900 illoo K 
{±111?-,° K) and 32 400±tJ!]!] K {±^„ K) at 68 (95) per cent 
confidence around the median. Nevertheless, it is encourag- 
ing that our procedure correctly distinguishes between mod- 
els with different temperatures. 

5.2 Measurement of To for SDSS J0818+1722 

Figure 7 displays the cumulative distribution for Tq, P{< 
TolDobs), derived for the line-of-sight SDSS J0818-hl722 at 
z = Q.We find a median temperature of To = 23 600±i§!5g K 
(igloo K), which lies in the middle of our simulation grid. 
The CPDF flattens towards high temperatures, indicating 
that temperatures of To > 33 000 K arc inconsistent with the 
observed Doppler parameter data. The lower limit on To is 
influenced by our choice of prior for To following Hi and 
He II reionisation, 13 500 K < To < 38 500 K {e.g. McQuinn 
et al. 2009). However, the gradient of the CPDF does begin 
to flatten towards the lower edge of the simulation grid, 
suggesting the probability of temperatures below 13 500 K 
is fairly low. 

5.3 Consistency check using proximity zone sizes 

In an analysis of observed proximity^ zone sizes, WBH08 
found the absolute sizes of the proximity zones in their 
RT models were systematically smaller than the observa- 
tional data. Based on their radiative transfer modelling of 
the proximity zones WBH08 further found that the size of 
the proximity zones should depend primarily on the back- 
ground H I photo- ionisation rate and the temperature of the 

^ WBH08 refer to these regions as highly ionised "near-zones", 
reflecting the fact that it can be diflicult to determine whether 
these transmission windows arc due to H II regions embedded in 
a significantly neutral IGM or a region of enhanced ionisation 
{i.e. the proximity effect) in an otherwise highly ionised IGM at 
2 > 6 (Bolton & Haehnelt 2007b; Maselli et ah 2007; Lidz et al. 
2007). However, in this work we mainly concentrate on the ionised 
IGM at z < 6, and so for simplicity we refer to these regions as 
proximity zones throughout. 
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Figure 7. Constraint on the temperature of the IGM at mean 
density at z = 6 in the proximity zone of the quasar SDSS 
J0818-I-1722. The solid curves show the cumulative distribution 
for To, P{< To I -Dobs), obtained using Eq. (6), while the shaded 
region indicates the 68 per cent confidence interval around the 
median. The numerical value in the top left of the panel corre- 
sponds to the recovered median temperature at mean density and 
the 68 (95) per cent confidence intervals. 



IGM, where the latter impacts on the ionisation state of 
hydrogen through the Hii recombination coefficient, ami oc 
T""'^. WBH08 therefore suggested that the IGM tempera- 
ture in quasar proximity zones at « ~ 6 may be as high as 
To ~ 40 000 K, perhaps resulting from the recent reionisa- 
tion of H I and He ii by a hard ionising spectrum. Imposing 
a high temperature of To ~ 40 000 K was found to resolve 
the discrepancy between the sizes of observed and simulated 
proximity zones. Our results, albeit for a single quasar prox- 
imity zone, are inconsistent with this interpretation. If our 
constraint is representative of the other z ~ 6 lines-of-sight, 
an alternative explanation would be required to explain the 
proximity zone sizes. 

In this paper we revisit the question of proximity zone 
size using our new suite of high resolution simulations. This 
analysis provides a consistency check of our temperature 
constraints. We have repeated the proximity zone size anal- 
ysis performed by WBH08 (see their section 5 for details) us- 
ing our model C simulations. This model has To — 21 500 K, 
similar to the median value derived for SDSS J0818-f 1722. 
The main difference between our work and WBH08 is 
the much higher resolution of our hydrodynamical simula- 
tions. WBH08 used a simulation with gas particle meisses of 
4.3 X lO'^h'^Mo and a box size of 60ft"^Mpc. This rather 
large volume and low resolution was chosen so that the 
largest halo in the simulation volume was reasonably mas- 
sive, with a mass of 2.7 x W^'^h~^MQ. Here we use much 
higher resolution hydrodynamical simulations which resolve 
the low density IGM (Bohon & Becker 2009). In this work 
we have demonstrated that mass resolution is more impor- 
tant than the halo mass for correctly resolving transmission 
from the Lya forest in the quasar proximity zone. We also 
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Figure 8. Proximity zone sizes (proper Mpc) measured from syn- 
thetic spectra (filled diamonds) compared to observational data 
taken from Carilli ct al. (2010). The simulated spectra are con- 
structed from model C, which has gas temperatures consistent 
with our constraint from the proximity zone of SDSS J0818-I-1722. 
Following Fan et al. (2006) and Carilli et al. (2010), the sizes have 
been rescaled to a common AB magnitude of M1450 = —27 by 
assuming the proximity zone size is proportional to N^/^. The 
inset displays the best fit linear relation to the observed and syn- 
thetic data, parameterised as -R = -Re — P(z — 6). The dotted 
and solid contours in the represent the 68, 84, 91 and 97 per cent 
bounds on single parameters for the observed and simulated data, 
respectively. 



now use the observational data compiled by Carilli et al. 
(2010), which include seven additional quasars and an im- 
proved determination of systemic redshifts. 

The results are displayed in Figure 8. The filled dia- 
monds correspond to the simulation data, while the open 
diamonds show the observational data compiled in Carilli 
et al. (2010). Following Fan et al. (2006) and CarilU et al. 
(2010), the sizes* have been rescaled to a common AB mag- 
nitude of Afi45o = —27 by assuming the proximity zone size 
is proportional to N^^'"'' , where N is defined by Eq. (1). Since 
our hydrodynamical simulation outputs are all at z = 6.01, 
we must also rescale the physical gas density in the simu- 
lations by (1 -I- z)^ to match the redshifts of the quasars, 
5.77 < z < 6.43. The inset shows the constraints on the 
amplitude and slope of the proximity zone evolution using 
the parameterisation R = Rs — p{z — 6), where Re is the 
value of the proximity zone size at « = 6 and /3 is the slope 
of the evolution. The dotted and solid contours represent 
the 68, 84, 91 and 97 per cent bounds on single parameters 
for the observed and simulated data, respectively. We find 
that the simulated data are in reasonable agreement with 
the observational data without having to arbitrarily impose 
high gas temperatures as in WBH08. This is because the 

Following the observational definition of Fan ct al. (2006), 
WBH08 define the extent of the proximity zones as the distance 
from the quasar at which the transmission falls below 10 per cent 
after smoothing by a top hat filter of width 20A. 



much higher mass resolution of our hydrodynamical simu- 
lation increases the transmission from the underdense re- 
gions. These dominate the transmission towards the edge of 
the proximity zone, increasing the extent of these regions. 
We have also verified that To > 33 000 K produces proxim- 
ity zone sizes which are too large by performing the same 
analysis using model G. Our temperature constraint from 
SDSS J0818-I-1722 is therefore consistent with the rapid size 
evolution of proximity zones at z ~ 6, providing a useful 
check of our modelling and results. As already noted by 
WBH08, the rapid evolution in the sizes of observed prox- 
imity zones also suggests that the IGM in the environment 
of bright quasars is already highly ionised by « ~ 6. 



6 IMPLICATIONS FOR THE REIONISATION 
HISTORY AROUND SDSS J0818+1722 

This work has focused on measuring the IGM temperature 
within close proximity to a quasar. Although this means we 
are unable to place constraints on the thermal state of the 
"average" IGM as previous studies have done at lower red- 
shift {e.g. Schaye et al. 2000; Ricotti et al. 2000; McDonald 
et al. 2001; Zaldarriaga 2002; Lidz et al. 2009), it does en- 
able us to study the IGM thermal history in the quasar's 
immediate vicinity (Miralda-Escude & Rces 1994; Theuns 
et al. 2002b; Mciksin et al. 2010). The proximity zone tem- 
perature will encode information on when the IGM near 
SDSS J0818-I-1722 was first reionised. The long timescale 
associated with adiabatic cooling in the low density IGM 
enables the gas to retain a "fossil record" of its initial tem- 
perature following reionisation (Hui & Gncdin 1997; Theuns 
et al. 2002a; Hui & Haiman 2003). Adopting an approach 
similar to that employed by Hui & Haiman (2003), we may 
therefore use our temperature measurement to place an up- 
per limit on the possible rcdshift of H i reionisation around 
SDSS J0818+1722. We pose the following question: "What 
is the maximum amount of time available for the IGM to 
cool following H i reionisation before the temperature be- 
comes inconsistent with our To constraint at 2 = 6?" 

To address this question wc shall compute the thermal 
evolution of a gas parcel at mean density in the presence 
of a power law ionising spectrum {e.g. Hui & Haiman 2003; 
Furlanetto & Oh 2009). Our code follows photo- ionisation 
and heating, collisional ionisation, radiative cooling, Comp- 
ton cooling and adiabatic cooling for six species (H i, H 11, 
He I, He 11, He iii, e~). We use the rates compiled in Bolton 
& Haehnelt (2007b) with the exception of the case B recom- 
bination and cooling rates of Hui & Gnedin (1997) and the 
photo-ionisation cross-sections of Verner et al. (1996). 

We consider two possibilities for the spectra of 
H I ionising sources (and hence IGM temperatures following 
H I reionisation) which should bracket the plausible range of 
possibilities. Following Hui & Haiman (2003), wc base the 
first model on the very massive, metal-free stellar spectrum 
presented by Bromm et al. (2001): Jv oc v for uai < v < 
I'Hei, oc for i/Hei ^ V < i^eii and Jv oc for 
V > i/Heii, where fi corresponds to the ionisation thresh- 
old for species i. This spectrum is very hard, and will re- 
sult in significant H i and He 11 photo-heating during reioni- 
sation. Note that since this model reionises He 11 at high red- 
shift, it is inconsistent with Hen reionisation completing at 
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Figure 9. Left: The IGM temperature at mean density, To, in a localised patch of the IGM as a function of redshift for two different 
instantaneous rcionisation models. The solid curve displays the temperature evolution following Hi and Hell reionisation at zjj = -^^Hc = 9.4 
by a very hard spectrum, representative of massive metal-free stars. The dashed curve assumes zn = 11.0 and corresponds to a much softer 
spectrum representative of population II sources. Note that Hell reionisation in this model is caused by the quasar itself at zjjc = 6.05 
(i.e. tq = lO'^ yr before z = 6). The filled square shows our constraint on the temperature in the proximity zone of SDSS J0818-I-1722. The 
thick (thin) error bars represent the 68 (95) per cent confidence interval around the median. Right: The cumulative probability distribution 
for the redshift of Hi reionisation around SDSS J0818+1722 we derive from the two different models. Our temperature measurement is 
consistent with Hi reionisation occurring at < 9.0 (zh < 11.0) at 68 (95) per cent confidence for the population II model, and zn < 8.4 
(zh < 9.4) for the metal-free spectrum. Note, however, that harder ionising spectra would increase both of these upper limits. 



lower redshift unless the spectrum softens significantly and 
Hem partially recombines again before z ~ 3. The second 
model, a more conventional population II spectrum, is based 
on the model of Leitherer ct al. (1999) for a galaxy of age 
500 Myr with a continuous star formation rate, a Salpeter 
IMF and metallicity Z = Q.2Zq. In this instance Jv oc 
below the Hi ionisation threshold. This soft spectrum does 
not significantly rcioiiisc He ii. Note that if the IGM is in- 
stead reionised by an earlier round of quasar activity, or if 
quasars are responsible for primarily driving H i reionisation 
(Volonteri & Gnedin 2009, but see Srbinovsky & Wyithe 
2007; Willott et al. 2009), the ionising spectrum and result- 
ing IGM thermal history will be intermediate between these 
two cases. 

To account for filtering of the ionising radiation as 
it propagates through optically thick gas, we assume the 
mean free path for ionising photons in a clumpy IGM has 
a frequency dependence oc v^'^ (Zuo & Phinney 1993; 
Miralda-Escude 2003). This hardens the intrinsic spectra 
by Os -> Os — 1.5, and is intermediate between the opti- 
cally thin case and the spectral modification expected for 
a uniform, optically thick IGM, <y.s — > — 3 (Abel & 
Haehnelt 1999). Lastly, we normalise both models to give 
an Hi photo-ionisation rate of Fhi = 1.9 x 10~^^ s~^, which 
is consistent with the upper limit inferred from the observed 
Lya forest opacity at 2: = 6 (Bolton & Haehnelt 2007c). 

Two example thermal histories are displayed in the 
left hand panel of Figure 9, along with our To constraint 
for the proximity zone of SDSS J0818-I-1722. Tlie tliick 
(thin) error bars correspond to the 68 (95) per cent con- 
fidence intervals around the median To. We assume instan- 
taneous H I and He 11 reionisation for the metal-free spec- 



trum at = ZYie = 9-4 (solid curve). The population II 
model (dashed curve) instead has 2h = 11.0. To account 

for Hen photo-heating by the hard spectrum of the quasar 
itself we also include He 11 rcionisation at zhc = 6.05 [i.e. 
tQ — 10^ yr before z — 6) for a quasar like spectrum with 
Ji/ oc !z~^'^ (hardened to Ji/ oc !z°). This leads to a tempera- 
ture boost during He 11 reionisation of ~ 8000 K, consistent 
with the results of our radiative transfer models. 

We stress that instantaneous reionisation is an unphysi- 
cal model for the ionisation and thermal history of the entire 
IGM during inhomogcneous Hi rcionisation. However, our 
temperature constraint is for the volume within 33 Mpc of a 
bright quasar. If jjliotohoating is the dominant heating pro- 
cess during the epoch of reionisation, we should nevertheless 
obtain an upper limit on the redshift when the bulk of the 
rcionisation occurred around the quasar in this approxima- 
tion. Numerical simulations indicate such biased regions arc 
amongst the first patches of the IGM to be reionised {e.g. 
Iliev et al. 2006; Shin et al. 2008; Finlator et al. 2009). A 
constraint on za within 33 Mpc SDSS J0818-I-1722 might 
then also give a constraint for the global onset of rcionisa- 
tion. However, a joint analysis with other observational con- 
straints including the Lya forest opacity (Fan et al. 2006) 
and the electron scattering optical depth (Dunkley et al. 
2009) will be required for a consistent interpretation of the 
temperature data (Choudhury & Ferrara 2006; Pritchard 
et al. 2009). In this work, our reionisation constraints apply 
to the IGM around SDSS JOS 18+ 1722 only. 

The right hand panel of Figure 9 displays the 
cumulative probability distribution for the redshift of 
H I reionisation, P{< zn), for these two models us- 
ing our z = 6 temperature measurement. Assuming 
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He II reionisation by the quasar boosts the IGM temperature 
by ~ 8000 K, we find the population II model is consistent 
with Hi reionisation occurring at z^i < 9.0 («h < 11.0) at 
68 (95) per cent confidence. However, a harder quasar spec- 
trum (qs < 1.5), and hence a larger temperature boost dur- 
ing He II reionisation, would weaken this constraint consider- 
ably. On the other hand, if H i and He ii reionisation around 
the quasar are driven primarily by a generation of metal free 
stars with hard spectra, we find 2h < 8.4 [zyi < 9.4) at 68 
(95) per cent confidence. 

We caution that the astrophysical uncertainties in the 
modelling presented here are significant; the exact spectrum 
of the ionising sources arc rather uncertain at 2: > 6. Harder 
ionising spectra would increase these upper limits. This 
study nevertheless demonstrates that even a single tempera- 
ture measurement at 2; = 6 can give an interesting constraint 
on the redshift of Hi reionisation around SDSS J0818-I-1722. 
Although computationally expensive, a larger set of sim- 
ulations will allow for a much more refined calibration of 
the measurement. Another obvious extension is to analyse 
additional lines-of-sight {e.g. Becker et al. 2005). Indepen- 
dent measurements along separate lines-of-sight will aid in 
reducing the statistical uncertainty on any averaged mea- 
surement. The amount of scatter in the temperature from 
one linc-of-sight to the next may also provide some insight 
into the topology and timing of H i reionisation globally. Im- 
proving existing constraints on the IGM thermal history at 
z < 4.5, as well as adding to constraints at higher redshift, 
will thus provide extremely valuable insight the epochs of 
both He II and H i reionisation. 



7 CONCLUSIONS 

In this paper we present the first direct measurements of 
the IGM temperature around a quasar &t z = We use a 
combination of high resolution hydrodynamical simulations 
combined with a line-of-sight radiative transfer implemen- 
tation to model the thermal state of the IGM in a quasax 
proximity zone. Previous theoretical studies have suggested 
that the IGM temperature close to a quasar will be sensi- 
tive to the prior ionisation state of hydrogen and helium, 
as well as the intrinsic quasar spectrum (Miralda-Escude & 
Rees 1994; Bolton & Haehnelt 2007b; Lidz et al. 2007). We 
demonstrate here that the Doppler parameter CPDF ob- 
tained from the quasar proximity zone provides a sensitive 
probe of the IGM temperature, enabling us probe the ther- 
mal and ionisation history of the IGM within the quasar's 
vicinity. 

Our observational data set consists of a single high 
resolution, high S/N Keck/HIRES spectrum of the quasar 
SDSS J0818-hl722 at z = 6. We perform identical Voigt pro- 
file analyses on the observed and synthetic data to obtain 
the Doppler parameter CPDF. The simulations, which use 
a range of self-consistent thermal histories, are used to cali- 
brate the line width measurements obtained from the obser- 
vational data. After verifying our method for obtaining tem- 
perature constraints, we proceed to obtain To = 23 600±s9[5[5 
K (±gi!]2 K) at 68 (95) per cent confidence around the me- 
dian for the proximity zone of SDSS J0818-I-1722. Due to 
the small size of the data set, the error bars are dominated 
by statistical uncertainty; we have verified that the most 



important systematic uncertainties remain small in compar- 
ison. We find our simulations are also in good agreement 
with the rapid evolution observed in proximity zone sizes 
at 2 ~ 6 when adopting temperatures consistent with our 
direct constraint. This provides a consistency check of our 
results, and further suggests that the IGM in the environ- 
ment of bright quasars is highly ionised by « ~ 6 (Wyithe 
et al. 2008; Maselli et al. 2009). 

Finally, we use our temperature constraint to explore 
the implications for the IGM reionisation history around 
SDSS J0818-M722 {e.g. Theuns et al. 2002a; Hui & Haiman 
2003). We consider two simple models for the thermal his- 
tory assuming that reionisation in the proximity zone oc- 
curred instantaneously. Assuming H i reionisation around 
SDSS J0818-I-1722 is driven by ionising sources with soft 
spectra, typical of population II stars, we infer zyi < 9.0 
(11.0) at 68 (95) per cent confidence if photo-heating is the 
dominant heating mechanism during the epoch of reionisa- 
tion. However, if the IGM is instead reioniscd by a popula- 
tion of massive metal-free stars, characterised by very hard 
ionising spectra, we obtain a tighter limit of zh < 8.4 (9.4). 
An H I reionisation epoch beginning at higher redshifts pro- 
duces temperatures which arc too low with respect to our 
constraint unless the H i ionising sources or the quasar itself 
have spectra significantly harder than typically assumed. 



ACKNOWLEDGEMENTS 

We thank Benedetta Ciardi and Peng Oh for valuable 
conversations during the course of this work. The hy- 
drodynamical simulations used in this work were per- 
formed using the Darwin Supercomputer of the Univer- 
sity of Cambridge High Performance Computing Service 
(http://www.hpc.cam.ac.uk/), provided by Dell Inc. using 
Strategic Research Infrastructure Funding from the Higher 
Education Funding Council for England. JSB acknowledges 
the support of an ARC Australian postdoctoral fellowship 
(DP0984947), and GDB thanks the KavU foundation for fi- 
nancial support. 

REFERENCES 

Abel, T. & Haehnelt, M. G. 1999, ApJ, 520, L13 
Alvarez, M. A. & Abel, T. 2007, MNRAS, 380, L30 
Becker, G. D., Ranch, M., & Sargent, W. L. W. 2007, ApJ, 

662, 72 

Becker, G. D., Sargent, W. L. W., & Ranch, M. 2005, in 
lAU CoUoq. 199: Probing Galaxies through Quasar Ab- 
sorption Lines, Williams, P., Shu, C.-G. & Menard, B. 
eds., p. 357 

Becker, R. H. et al. 2001, AJ, 122, 2850 

Bolton, J., Moiksin, A., & White, M. 2004, MNRAS, 348, 
L43 

Bolton, J. S. & Becker, G. D. 2009, MNRAS, 398, L26 
Bolton, J. S. & Haehnelt, M. G. 2007a, MNRAS, 381, L35 
Bolton, J. S. & Haehnelt, M. G. 2007b, MNRAS, 374, 493 
Bolton, J. S. & Haehnelt, M. G. 2007c, MNRAS, 382, 325 
Bolton, J. S., Oh, S. P., & Furlanetto, S. R. 2009, MNRAS, 
395, 736 



© 0000 RAS, MNRAS 000, 000-000 



The IGM temperature around a quasar at z — 6 15 



Bromm, V., Kudritzki, R. P., & Loeb, A. 2001, ApJ, 552, 
464 

Carilli, C.L. et al. 2010, ApJ, in press (arXiv:1003.0016) 
Cen, R., McDonald, P., Trac, H., & Loeb, A. 2009, ApJ, 

706, L164 

Choudhury, T. R. & Ferrara, A. 2006, MNRAS, 371, L55 

Croton, D. J. 2009, MNRAS, 394, 1109 

Dijkstra, M., Lidz, A., & Hui, L. 2004, ApJ, 605, 7 

Dunkley, J, et al. 2009, ApJS, 180, 306 

Eisenstein, D. J. & Hu, W. 1999, ApJ, 511, 5 

Fan, X., Narayanan, V. K., Strauss, M. A., White, R. L., 

Becker, R. H., Pcntcricci, L., & Rix, H.-W. 2002, AJ, 123, 

1247 

Fan, X. et al. 2006, AJ, 132, 117 

Faucher-Giguere, C.-A., Prochaska, J. X., Lidz, A., Hern- 

quist, L., & Zaldarriaga, M. 2008, ApJ, 681, 831 
Finlator, K., Ozel, F., Dave, R., &: Oppenheimer, B. D. 

2009, MNRAS, 400, 1049 
Furlanetto, S. R. & Oh, S. P. 2008, ApJ, 681, 1 
Furlanetto, S. R. & Oh, S. P. 2009, ApJ, 701, 94 
Gunn, J. E. & Peterson, B. A. 1965, ApJ, 142, 1633 
Haardt, F. & Madau, P. 2001, in Glusters of Galaxies 
and the High Redshift Universe Observed in X-rays, Neu- 
mann, D. M. & Tran, J. T. V. ed., preprint (arXiv:astro- 
ph/0106018) 

Hachnclt, M. G., Natarajan, P., & Rees, M. J. 1998, MN- 
RAS, 300, 817 

Haehnelt, M. G. & Steinmetz, M. 1998, MNRAS, 298, L21 

Hui, L. & Gnedin, N. Y. 1997, MNRAS, 292, 27 

Hui, L. & Haiman, Z. 2003, ApJ, 596, 9 

Hui, L. & Rutledge, R. E. 1999, ApJ, 517, 541 

Ilicv, I. T., Mellema, G., Pen, U.-L., Mcrz, H., Shapiro, 

P. R., & Alvarez, M. A. 2006, MNRAS, 369, 1625 
Kelson, D. D. 2003, PASP, 115, 688 

Kim, T.-S., Bolton, J. S., Viel, M., Haehnelt, M. G., & 
Carswell, R. F. 2007, MNRAS, 382, 1657 

Kirkman, D. ct al. 2005, MNRAS, 360, 1373 
Komatsu, E. et al. 2009, ApJS, 180, 330 
Leitherer, C. et al. 1999, ApJS, 123, 3 
Li, Y. ct al. 2007, ApJ, 665, 187 

Lidz, A., Fauchcr-Gigucrc, C., Dall'Aglio, A., McQuinn, 
M., Fechner, C, Zaldarriaga, M., Hcrnquist, L., & Dutta, 
S. 2009, ApJ, submitted (arXiv:0909.5210) 

Lidz, A., McQuinn, M., Zaldarriaga, M., Hernquist, L., & 
Dutta, S. 2007, ApJ, 670, 39 

Madau, P. & Meiksin, A. 1994, ApJ, 433, L53 

Martini, P. 2004, in Coevolution of Black Holes and Galax- 
ies, Ho, L. C. ed., Cambridge Univ. Press, Cambridge, 
p.169 

Maselli, A., Ferrara, A., & Gallerani, S. 2009, MNRAS, 
395, 1925 

Maselli, A., Gallerani, S., Ferrara, A., & Choudhury, T. R. 

2007, MNRAS, 376, L34 
McDonald, P., Miralda-Escudc, J., Raucli, M., Sargent, 

W. L. W., Barlow, T. A., & Ccn, R. 2001, ApJ, 562, 52 
McDonald, P., Miralda-Escude, J., Ranch, M., Sargent, 

W. L. W., Barlow, T. A., Cen, R., & Ostriker, J. P. 2000, 

ApJ, 543, 1 

McQuinn, M., Lidz, A., Zaldarriaga, M., Hernquist, L., 
Hopkins, P. F., Dutta, S., & Faucher-Giguere, C.-A. 2009, 
ApJ, 694, 842 



Meiksin, A., Tittley, E. R., & Brown, C. K. 2010, MNRAS, 
401, 77 

Mesinger, A. & Haiman, Z. 2007, ApJ, 660, 923 
Miralda-Escude, J. 2003, ApJ, 597, 66 
Miralda-Escude, J. & Rees, M. J. 1994, MNRAS, 266, 343 
Olive, K. A. & Skillman, E. D. 2004, ApJ, 617, 29 
PawUk, A. H., Schaye, J., & van Scherpenzeel, E. 2009, 

MNRAS, 394, 1812 
Peeples, M. S., Weinberg, D. H., Dave, R., Fardal, M. A., 
& Katz, N. 2009, MNRAS, submitted (arXiv:0910.0256) 
Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flan- 
nery, B. P. 1992, Numerical recipes in FORTRAN. The art 
of scientific computing, ed. W. H. Press, S. A. Teukolsky, 
W. T. Vetterling, & B. P. Flanncry 
Pritchard, J. R., Loeb, A., & Wyithe, J. S. B. 2009, MN- 
RAS, submitted (arXiv:0908.3891) 
Ranch, M. et al. 1997, ApJ, 489, 7 
Rcichardt, C. L. ct al. 2009, ApJ, 694, 1200 
Ricotti, M., Gnedin, N. Y., & ShuU, J. M. 2000, ApJ, 534, 
41 

Schaye, J., Theuns, T., Leonard, A., & Efstathiou, G. 1999, 
MNRAS, 310, 57 

Schaye, J., Theuns, T., Ranch, M., Efstathiou, G., & Sar- 
gent, W. L. W. 2000, MNRAS, 318, 817 

Shin, M., Trac, H., & Cen, R. 2008, ApJ, 681, 756 

Songaila, A. 2004, AJ, 127, 2598 

Springcl, V. 2005, MNRAS, 364, 1105 

Springcl, V. & Hcrnquist, L. 2003, MNRAS, 339, 289 

Srbinovsky, J. A. & Wyithe, J. S. B. 2007, MNRAS, 374, 
627 

Teller, R. C, Zheng, W., Kriss, G. A., & Davidsen, A. F. 

2002, ApJ, 565, 773 

Theuns, T., Leonard, A., Efstathiou, G., Pearce, F. R., & 

Thomas, P. A. 1998, MNRAS, 301, 478 
Theuns, T., Schaye, J., & Haehnelt, M. G. 2000, MNRAS, 

315, 600 

Theuns, T., Schaye, J., Zaroubi, S., Kim, T., Tzanavaris, 
P., & Carswell, B. 2002a, ApJ, 567, L103 

Theuns, T., Zaroubi, S., Kim, T.-S., Tzanavaris, P., & Car- 
swell, R. F. 2002b, MNRAS, 332, 367 

Trac, H., Cen, R., & Loeb, A. 2008, ApJ, 689, L81 

Tytlcr, D. et al. 2004, ApJ, 617, 1 

Valagcas, P., Schacffcr, R., & Silk, J. 2002, A&A, 388, 741 
Verner, D. A., Ferland, G. J., Korista, K. T., & Yakovlev, 

D. G. 1996, ApJ, 465, 487 
Volonteri, M. & Gnedin, N. Y. 2009, ApJ, 703, 2113 
Volonteri, M. & Rees, M. J. 2006, ApJ, 650, 669 
Walter, F., Carilli, C, Bcrtoldi, F., Menten, K., Cox, P., 

Lo, K. Y., Fan, X., & Strauss, M. A. 2004, ApJ, 615, L17 
White, R. L., Becker, R. H., Fan, X., & Strauss, M. A. 

2003, AJ, 126, 1 

Willott, C. J. et al. 2009, AJ, in press (arXiv:0912.0281) 
Wyithe, J. S. B., Bolton, J. S., & Haehnelt, M. G. 2008, 

MNRAS, 383, 691 
Zaldarriaga, M. 2002, ApJ, 564, 153 
Zuo, L. & Phinney, E. S. 1993, ApJ, 418, 28 



© 0000 RAS, MNRAS 000, 000-000 



